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Abstract. The nonlinear integral equation P(x) = j dyw(y) P(y) P(x + y) is 
investigated. It is shown that for a given function w(x) the equation admits an infinite 
set of polynomial solutions P(x). For polynomial solutions, this nonlinear integral 
equation reduces to a finite set of coupled linear algebraic equations for the coefficients 
of the polynomials. Interestingly, the set of polynomial solutions is orthogonal with 
respect to the measure xw(x). The nonlinear integral equation can be used to specify 
all orthogonal polynomials in a simple and compact way. This integral equation 
provides a natural vehicle for extending the theory of orthogonal polynomials into 
the complex domain. Generalizations of the integral equation are discussed. 



There are many ways to specify uniquely a set of orthogonal polynomials. One can 
specify the domain (a, (3) and the measure with respect to which the polynomials are 
orthogonal and then use the cumbersome Gramm-Schmidt orthogonalization procedure 
to construct the polynomials. For example, for the domain (—1, 1) and measure 
(1 — x 2 ) -1 / 2 , the Gramm-Schmidt procedure yields the Chebyshev polynomials T n {x). 
Alternatively, one can specify a recursion relation. The linear recursion relation 
T n+ i(x) = 2xT n (x) — T n _i(x) along with the initial conditions T Q (x) = 1 and Ti(x) = x 
again produces the Chebyshev polynomials. Another approach is to give the differential- 
equation eigenvalue problem satisfied by the polynomials. The Chebyshev polynomials 
T n (x) satisfy the eigenvalue equation (1 — x 2 )y"(x) — xy'(x) + n 2 y(x) = 0. Stating 
the generating function is yet another way to specify a set of polynomials. For the 
Chebyshev polynomials the generating function 



uniquely defines these polynomials. 

In this paper we propose a simple and compact way to specify a set of orthogonal 
polynomials in terms of a nonlinear integral equation. Consider the nonlinear integral 
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equation 

P(x)= / dyw(y)P(y)P(x + y). (1) 

J a 

The integration limits a and (3 as well as the function w(x) are arbitrary except for the 
restriction that w(x) must have a nonvanishing integral, f dxw(x) ^ 0. Therefore, 
without loss of generality, we may assume that w is normalized: 

P 

dxw(x) = 1. 

Note that if we choose a = — oo, (3 = oo, and w(y) = e iy , ([I]) reduces to the heavily 
studied equation for the Wigner function [H[2]. The original motivation for considering 
this nonlinear integral equation is that it describes a stochastic process in which two 
random variables are subtracted to create a new one. The steady-state probability 
distribution for this random variable satisfies the integral equation (CQ) with a = 0, 
(3 = oo, and w (y) = 2. 

Even though this integral equation is nonlinear, its polynomial solutions can be 
found analytically. By inspection one can see that there is the trivial solution Pq(x) = 1. 
However, there are also infinitely many other polynomial solutions because ([1]) has two 
remarkable properties. 

First, if we seek a solution in the form of an arbitrary polynomial of degree n, 

n 

P n (x) = ^a nik x k , (2) 

fc=0 

the nonlinear integral (OQ) preserves the degree of the polynomial. For example, if we 
substitute an arbitrary linear polynomial P\{x) = a 10 + ai ix or an arbitrary quadratic 
polynomial i^^) = ^2,o + a 2,i% + a2,2% 2 into ([T|), we obtain 

ai,o + a>i,i% = (Pi(y) [ai,o + ai,iV + 

and 

«2,o + d2,ix + a 2) 2X 2 = (P-ziv) [«2,o + 02,12/ + a 2,ix + a 2 ^y 2 + 2a 2 , 2 xy + a 2 , 2 y 2 ] ), 
where we have introduced the notation 

(/)= I dxw(x)f(x). (3) 



In general, for an arbitrary polynomial of degree n both the left and right sides of (Op) 
are polynomials of the same degree n. 

Substituting a polynomial of degree n into the right side of (pQ) and performing the 
integration, we obtain a polynomial of degree n. 

Second, for the case of polynomial solutions, the nonlinear equation <^ reduces to 
a system of coupled linear equations, which are obtained by equating like powers of x. 
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For the case n = 1 above, equating the coefficient of x 1 on both sides of the equation 
gives a lt i = a 11 (P 1 (?/)). We require that the polynomial have degree one, a^i ^ 0, so 

(P 1 (x)) = 1, (4) 

where we have replaced y by x. Next, we equate the coefficients of x° and get 
a i,o — a i,o{Pi(y)) + a i,i{y Pi (y))- Substituting (jlj) into this equation, we obtain 

(xP 1 (x))=0. (5) 

Thus, the nonlinear equation ([1]) reduces to a system of two coupled inhomogeneous 
linear equations for the two unknowns ai j0 and a^i. The integral equation reduces to a 
linear algebraic system only for polynomials because, as we see in the above calculation, 
the derivation relies on the fact that the power series in x truncates at a finite order. 

Let us repeat this procedure for quadratic polynomials: Starting at the highest 
power x 2 , we find that 02,2 = 02,2(^2(2/)) and hence we recover with the subscript 1 
replaced by 2. Setting the coefficients of x 1 equal gives a 2i i = ci2, 1(^2(2/)) + 2a 2 ,2 (2/^2(2/)) 
and we recover ([5]) with the subscript 1 replaced by 2. The coefficient of x° gives 
02,0 = CL 2 ,o(p2(y)) + 02,1(2/^2(2/)) + CL 2 ,2(y 2 P2(y)} ■ Therefore, there is now the third 
equation 

(x 2 P 2 {x)) =0. (6) 

The general pattern is clear. An nth-degree polynomial is a solution of the integral 
equation if and only if the following set of n + 1 linear equations is satisfied: 

(x k P n (x)) = 5 kfi (k = 0,l,...,n). (7) 

The linear equations ([7]) imply that the polynomials P n {x) are orthogonal with 
respect to the measure 

g(x) = xw(x). (8) 

If P m (x) = YHk=o a m,kX k is a polynomial of degree m < n, then from ([7]) we conclude 
that 

m m+1 
(xP n P m ) = J2 a ™,k(x k+1 Pn(x)} = J2 a ^-l(x k Pn(x)} = 0. (9) 
fc=0 k=l 

This is the main result of our paper. We have shown that the nonlinear integral equation 
(Jl]) admits an infinite set of polynomial solutions. This set is unique; there is one and 
only one polynomial solution of degree n. These polynomials are orthogonal. 

We comment that this integral formulation of orthogonal polynomials is general. 
The integration limits and the function w(x) are arbitrary apart from the restriction 
that w(x) have a positive integral J dxw(x) > 0. 

We mention three classical orthogonal polynomials |3j that can be generated using 
this approach: 

• Generalized Laguerre polynomials Ln\x) using a = 0, (3 = oo, and w(x) = 
x 1 ~ 1 e~ x for all 7 > 1; 
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\P-Q 



• Jacobi polynomials G n (p,q,x) using a = 0, (3 = 1, and w(x) = x q ~ 2 (l 
for q > 1 and p — q > — 1; 

• Shifted Chebyshev polynomials of the second kind Z7* (x) using a = 0, (3 = 1, and 

= (1 -x) 1 /^" 1 / 2 . 

The equations (JTj) can be written compactly as n + 1 simultaneous linear equations 
for the n + 1 coefficients a r , 



(10) 



a nJ m k+j = 4,o (fe = 0, 1, . . . , n). 

j=0 

Here m n = (x n ) are the moments of w(x) and by assumption m = 1. Thus, using 
Cramer's rule we can express the polynomial solutions explicitly as ratios of determinants 
[I]. We define the matrices 



( 1 x 

mi m 2 

\ m n m m+1 



x n \ 
m n+1 

m 2n J 



( m mi 

mi m-i 

\ m n m m+l 



m n \ 
m n+ i 

m 2n J 



so that B n = (Ap). The polynomials P n (x) are then given by 

det A n 



PJx) 



det B r 



The normalization of the polynomials can be given in the form 

(x P n {x) P m {x)) 0~n,m G n - 

The normalization factors G n are ratios of determinants of matrices 
(detC n )(detC n+1 ) 



G n 



(detB n ^ 2 



where the matrix C n is given by 



Cn 



( mi m 2 

m 2 W3 



m n +i 
m 2n J 



\ m n m m+ i ■ 

For completeness, we give the first three polynomials explicitly: 
P (x) = 1, 



Pi(x) 

P2(x) 



m 2 — xmi 
m 2 — m\ ' 



(m 2 m l i — m^) + {m 2 m^—mim4)x+ (mim^ — m\)x 2 
m i (m 2 — m\) — m\ + 2m 1 m 2 m 3 — m 2 



(11) 
(12) 

(13) 
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Also, we give the corresponding normalization factors G n : 
G = m 1 , 
Gi = mi 



1\2 ' 



(777,2 — m 



[m4(m2 — mf) — Trig + 2mi7Tl27Tl3 — 7772 12 



I ■ 1 I 1 ./ill. 1 > 1 — ^ 1 > ! ~ ' ) t . I III III i 1 ( 1 

Note that f|TT| - [T3|) are only valid if det B n 7^ 0. This condition holds when the measure 
g(x) = xw(x) is positive 01a a < x < j3 [S1E]- 

The integral equation ([1]) specifies all possible orthogonal polynomials. When the 
weight function g(x) with respect to which the polynomials are orthogonal does not 
vanish at x = and consequently w(x) = g(x)/x is singular at x — 0, the path of 
integration from a to (3 may be taken in the complex plane to avoid the singularity at 
the origin. Using a complex integration path, all steps leading to ([9]) are valid. 

For example, consider the Legendre polynomials for which a = —1, (3 = 1, and g(x) 
is a constant. There are an infinite number of topologically distinct integration paths 
that connect —1 to 1. These paths are characterized by their winding numbers. For 
definiteness, we choose a path that goes from —1 to 1 in the positive (counterclockwise) 
direction and does not encircle the origin. On this path J dx/x = in and hence to 
maintain the normalization J^dxw(x) = 1 we use w(x) = l/(i7rx). The moments 
m n = (x n ) are mo = 1, ttii = 2/ (in), 1712 = 0, 7773 = 2/(3z7r), .... From the moment 
formulas ffTTl) . we obtain 

P (x) = 1, 

VK 

P x [x) = —x, 
p 2 (x) = 1 - 3x 2 , 

P 3 (x) = ^(3x-5x 3 ), (14) 

o 

and so on. These polynomials are the standard Legendre polynomials, except that 
the odd polynomials have an imaginary multiplicative factor that is determined by the 
winding number of the integration path. Of course, these polynomials are solutions of 
the linear equations ([7j). 

To summarize, while the usual theory of orthogonal polynomials is formulated in 
terms of real integrals, the integral equation (CD) provides a simple and natural framework 
to extend the theory of orthogonal polynomials into the complex domain. In doing 
so we discover an interesting connection between the polynomial coefficients and the 
topological winding number of the integration path. 

There are many ways to generalize the integral equation (CD): 
1. Multiplicative argument. If we replace the term P(x + y) in the original nonlinear 
integral equation ([1]) by P(xy), we obtain a new class of nonlinear equations: 

rP 

P(x)= / dyw(y)P(y)P(xy). (15) 
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Each of these nonlinear integral equations also has an infinite number of polynomial 
solutions. 

Again, there is the constant solution Pq(x) = 1. To find other solutions we 
substitute a polynomial of degree n and then equate coefficients of x k on the left and 
right sides, starting with k = n. The nonlinear integral equation reduces to a set of 
n + 1 inhomogeneous equations of the form 

a U)k {x k P n {x)) = a nM (k = 0,1, . . . ,n). (16) 

However, unlike the previous case, the equations are quadratic and there are now 
2™ _1 solutions because each of the coefficients a n ^ can be either zero or nonzero for 
k — 0, . . . , n — 1. For each nonzero coefficient the linear equation 

(x k P n {x)) = 1 (17) 

holds for all k for which a n j, ^ 0. 

In one special class of solutions all of the coefficients are nonzero and the 
polynomials are orthogonal with respect to the measure 

g(x) — (1 — x) w(x). (18) 

To verify this, we take a polynomial of degree m < n and observe that 

m 

((1 - x)P n P m ) = ^a m , fc {(x k P n ) - (x k+1 P n )) = 0. 

fc=0 

This equation is valid because all of the terms in the parentheses vanish by virtue of 
(flTl) . This measure is relevant for the class of polynomials for which < a < (3 < 1 . 
Here are two classical orthogonal polynomials that can be generated in this way: 

• Jacobi polynomials G n (p,q,x) using a = 0, (3 = 1, and w(x) = (1 — x) p ~' ? ~ 1 £' 3 ~ 1 
with p — q > and q > 0; 

• Shifted Chebyshev polynomials of the second kind U*(x) using a = 0, (3 = 1, and 
w [x) = x 1 ' 2 (l-x)- l l 2 . 

In another class of solutions, the coefficients alternate between zero and nonzero so 
that the polynomials alternate between definite even and odd parity. The polynomials 
are orthogonal with respect to a measure whose moments 

fJln = I( X n_ X n+2)[l + (_l)n]. (19) 

Thus, all the even moments of the measure are positive and all the odd moments vanish. 
It is easy to show that polynomials of similar parity are orthogonal with respect to the 
measure (1 — x 2 )w(x). Polynomials of dissimilar parity are orthogonal because their 
product is an odd polynomial, and the odd moments vanish (TT9]) . 
2. Linear arguments. For the integral equation with a linearly shifted argument, 

rP 

P(x)= / dyw(y)P(y)P(x + a + by), (20) 

J a 

where b ^ is an arbitrary constant, polynomials of degree n are solutions when 
((a + bx) k P n (x)) = 5kfi- This implies that the polynomials are orthogonal with respect 
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to the measure g(x) = (a + bx)w(x). Note that when a = 0, the polynomials are 
identical to those generated by the original integral formula ([I]). Curiously, for a = 1 
and b = — 1 we recover the polynomials generated by the integral equation ( jl~5l) . 
3. Functional arguments. The integral equation 



where f(x) is a nonconstant function, has polynomial solutions of degree n when 
([f( x )] k Pn{x)) = 5k,o- These polynomials are not necessarily an orthogonal set. 
Nevertheless, the polynomial P n (x) is orthogonal to the function P m [f(x)\ when m < n 
with respect to the measure g(x) = f(x)w(x). 
4. Arbitrary functions. The integral equation 



where f(x) is a nonconstant function, has polynomial solutions of degree n when 
{x k f[P n (x)]) = 5fe )0 . This implies that the function f[P n (x)] is orthogonal to the 
polynomial P m {x) with respect to the measure g(x) = xw(x). 

5. Further generalizations. It is worth considering what happens when the function 
w(x) is singular on the interval a < x < (3. Also, an obvious way to generalize ([1]) is to 
iterate it, and thereby to obtain integral equations that are cubic, quartic, and so on. 
Furthermore, one can study the properties of nonpolynomial solutions to ([IT) . We have 
found many such solutions. Finally, one can generalize (0Q) to multidimensional integrals 
and study the properties of the resulting multivariate polynomial solutions. 

In summary, we have shown that all sets of orthogonal polynomials are solutions of 
nonlinear integral equations. For polynomial solutions these nonlinear equations reduce 
to simultaneous linear equations for the coefficients of the polynomials. The measure 
with respect to which the polynomials are orthogonal depends on the form of the integral 
equation and on the integration measure. The nonlinear integral equations discussed 
here provide a simple and compact way to define a set of orthogonal polynomials and 
also provide a framework for extending the general theory of orthogonal polynomials 
into the complex domain. 
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